Memory beyond memory in heart beating: an efficient way to detect pathological 

conditions 



O 

o 



O 
CO 



o 



CO 



I 

C 
O 
O 



> 

in 

o 



a 



i 

c 
o 
o 



X 
S3 



P. Allegrini 1 , P. Grigolini 2 ' 3 ' 4 , P. Hamilton 5 , L. Palatella 3 , G. Raffaelli 3 

1 Istituto di Linguistica Computazionale del Consiglio Nazionale delle Ricerche, 

Area della Ricerca di Pisa-S. Cataldo, Via Moruzzi 1, 56124, Ghezzano-Pisa, 

2 Center for Nonlinear Science, University of North Texas, P.O. Box 311427, Denton, Texas, 76203-1427 

3 Dipartimento di Fisica dell'Universita di Pisa and INFM Piazza Torricelli 2, 56127 Pisa, Italy 

4 Istituto di Biofisica del Consiglio Nazionale delle Ricerche, 

Area della Ricerca di Pisa-S. Cataldo, Via Moruzzi 1, 56124, Ghezzano-Pisa, Italy 

5 Center for Nonlinear Science, Texas Woman's University, P.O. Box 425498, Denton, Texas 76204 

(February 1, 2008) 



We study the long-range correlations of heartbeat fluctua- 
tions with the method of diffusion entropy. We show that this 
method of analysis yields a scaling parameter 8 that appar- 
ently conflicts with the direct evaluation of the distribution of 
times of sojourn in states with a given heartbeat frequency. 
The strength of the memory responsible for this discrepancy 
is given by a parameter e 2 , which is derived from real data. 
The distribution of patients in the (8, e )-plane yields a neat 
separation of the healthy from the congestive heart failure 
subjects. 

PACS numbers: 87.19.Hh, 05.45.Tp, 05.40. Fb 

The analysis of time series of physiological significance 
is currently done using the paradigm of anomalous scal- 
ing JL|. This letter, resting on this paradigm, aims at 
showing that the entropy of a diffusion process generated 
by a physiological time series according to the prescrip- 
tions of Rcfs. y||] yields a scaling exponent that depends 
only on genuine events, namely, events whose occurrence 
time is unpredictable. We call this method of analysis 
Diffusion Entropy (DE) method and, by means of simple 
dynamic models, we prove it to be insensitive to pseudo 
events, namely, to events whose occurrence times are cor- 
related to those of earlier events. 

Let us consider first a dynamic model that generates 
events (really random events). This is given by 



x = $(#) > 0, 



(1) 



where x denotes the coordinate of a particle, moving 
within the interval / = [0, 1], from the left to the right, 
with times of arrival at x = 1 determined by Eq. (|]J) and 
by the initial condition. When the particle reaches the 
right border of /, it is injected back to a new initial con- 
dition selected with uniform probability on /. Conse- 
quently, the times of arrival at x = 1, t\, ..,tj.., are a fair 
example of real events. It is straightforward to prove that 
the choice <&(x) — nx z , with z > 1 and k > 0, yields for 
the waiting times Tj = ti — t%-i the following distribution 
density 



with fj, = z/(z — 1) and T = (/i — 1)/k. Note that the 
mean waiting time (r) is determined by T through (r) = 
T/(/i — 2). Let us convert now the time series {tj} into 
a random walk. We select the rule [g| that makes the 
random walker move, always in the same direction and by 
a step of constant intensity, only when an event occurs. 
This means that the sequence {r^} is converted into a 
sequence of 0's and l's as follows. We create a sequence 
of patches, one patch for each n, with a width given 
by the integer part of r^; we fill the sites of each patch 
with 0's and we signal the border between two nearest- 
neighbor patches with l's. Then, the resulting sequence 
is converted into many trajectories of a given length / 
with the method of the moving window. A window of size 
I moves along the sequence and for any window position, 
the portion of the whole sequence spanned by the window 
is regarded as a single trajectory of length I. All these 
trajectories are assumed to start from the origin, and 
are used to create a diffusion distribution, at time I. If 
there is scaling, with scaling parameter 8, the entropy 
S(l) takes the form @,| 



S(Ji)=A + 6bx(l). 



(3) 



Consequently, the rate of entropy increase in logarithmic 
time scale is the value of the scaling parameter we are 
looking for. The DE method has been widely discussed 
in earlier papers p,H and it has been pointed out that, 
without using any form of detrending, it yields the cor- 
rect scaling parameter 8. In some cases jlj this method 
detects the Levy scaling that would imply an infinite vari- 
ance and, consequently, would be incompatible with the 
adoption of methods based on variance. In the specific 
case of the model of Eq. dh) the adoption of this method 
of analysis yields a scaling parameter 5 that, in the case 
where 2 < /x < 3, fits the theoretical prediction 



6 = 



1 



M-l' 



(4) 



^(r) = (/i - 1) 



rpn- 



{T + r y 



(2) 



while for /i > 3 one gets 8 = 0.5. This is a proof that the 
direct evaluation of the power law index /i is equivalent to 
detecting the scaling 8 by means of the DE method. It is 
well known B that a Markov master equation, namely a 



stochastic process without memory, is characterized by a 
waiting time distribution ip(r) with an exponential form, 
thereby implying that a marked deviation from the expo- 
nential condition is a signature of the presence of mem- 
ory. This is the kind of memory that is usually associated 
with the detection of an anomalous scaling parameter, 
namely 8 ^ 0.5. We refer ourselves to this memory as 
memory of the first type. To prepare the ground for the 
kind of memory that is the focus of this paper, we have 
to discuss a dynamic model generating both events and 
pseudo events. For this purpose let us consider a two- 
variables model. The first equation, referring to the first 
variable, is given by Eq.([y), and the second equation, 
concerning the new variable y, is given by 



y = x(y) > o. 



(5) 



The variables x and y are the coordinates of two par- 
ticles, both moving in the interval /, always from the 
left to the right. The initial conditions of the variable 
y are always chosen randomly. The initial conditions of 
x, on the contrary, are not always chosen randomly, but 
they are only when the variable y reaches the border at 
least once, during the sojourn of x within the interval. 
Let us consider the sojourn time interval [tj,fi+i]. If in 
this time interval the variable y remains within the in- 
terval, without touching the right border, then we set 
x(ti + i) — x(U). This is a pseudo event. Thus, the se- 
quence {ti} is a mixture of events and pseudo events. Let 
us consider the case where x(y) = k'y z with z' > 1 and 
k' > 0, so as to produce the power index // = z 1 j(z! — 1), 
with /i' > 2, a property of real events. Let us set the 
condition (t) x <C (t) v . In this case it is straightforward 
to show that the waiting time distribution of x is still 
given by Eq.(EI). However, the power index [i is not more 
a reflection of real events. Fig.^j reveals a very attractive 
property: the DE now yields 8 = l/(// — 1), which is 
very different from the theoretical prediction of Eq.(|J) 
that makes /i > 3 yield 8 — 0.5. Note that in the case 
of Fig.0 fi = 5 and // = 2.41. The breakdown of Eq.(f|) 
is referred to by us as memory beyond memory effect. 
In fact, the existence of pseudo events implies correla- 
tion among different times of the series {tj}, and thus a 
memory of earlier events. This memory is annihilated by 
shuffling the order of these times, as we do with a bunch 
of cards. In fact, as shown by the inset of Fig.pl we see 
that shuffling has the effect of yielding 8 = 0.5, in accor- 
dance with the theoretical prescription of Eq.(Q). It is 
impressive that the scaling detected by the DE method 
does not depend on the pseudo events, but only on the 
hidden events that would be invisible to an analysis based 
on the direct evaluation of ip{ T )- 

Let us apply this model to the real data taken from 
pL We apply our technique to 33 long-time ECG 
records (about 20 hours each), 18 healthy and 15 showing 
pathological behavior related to Congestive Heart Fail- 
ure (c.h.f.). Following Ref. 0, we refer to all the ECG 



records of the MIT-BIH Normal Sinus Rhythm Database 
and of the BIDMC Congestive Heart Failure Database, 
for the healthy and the c.h.f. patients, respectively. 

The data under study are time series of the kind of 
that illustrated in Fig. 2, where the length of the vertical 
lines expresses T(i) = ti— U-i as a function of the integer 
number i. The integer number i denotes the i-th heart 
beating of an electrocardiogram, and ti is the time at 
which the R wave of this heart beat occurs. 
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FIG. 1. DE for two-variables model as a function of time. 
The squares correspond to fc'=0.018, z'=1.83 while the dia- 
monds to fc'=0.011, z'—1.71. For both curves k-OA, 2=1.25. 
In the inset: the same curves after shuffling, the straight line 
slope is 0.5. 



We make these data suitable for the illustration of the 
memory beyond memory effect as follows. We adopt a 
procedure illustrated with the help of Fig.g. The vertical 
axis, concerning the variable T(i), is divided into many 
cells of a given size AT. Thus the (T(i),z)-plane is di- 
vided into many horizontal strips with a constant width 
equal to AT. This coarse-graining prescription yields 
the thick line of Fig|2| The curve corresponds to many 
horizontal intervals separated by vertical up and down 
jumps. The widths of these horizontal intervals define a 
sequence of numbers Ti that is the object of our statistical 
analysis. To make this analysis as efficient as possible we 
have to make a proper choice of the value of AT, since an 
excessively small value would produce too many pseudo 
events and an excessively large would yield poor statis- 
tics. The results of our statistical analysis were proven 
to be insensitive to changing AT over a wide range of 
values. We therefore assign to AT the mean value of this 
range, which turns out to be AT = 1/30 sec. 

The events under study refer to the jumps from one 
to another strip. To assess whether they are events or 
pseudo events, we have to compare the waiting time dis- 
tribution tp(r) to the scaling detected by means of the 
DE method. For the sake of statistical accuracy we de- 
cided to evaluate the probability of finding waiting times 



larger than a given value r. This is the function ^(r) 
defined by ^(r) = /. dT0(r). The results illustrated in 
the inset of Fig|| imply the Brownian scaling S = 0.5. 
In fact, the function \&(r) of the heart failure subject is 
a stretched exponential and the healthy subject yields 
\x = 3.9. The same Brownian condition applies to all the 
subjects. On the contrary the DE method yields for the 
healthy subjects the mean value S — 0.82 ± 0.04 and for 
the heart failure subject 6 — 0.71 ± 0.06. It is interesting 
to notice that Fig.0 refers to the same subjects as those 
of the inset of FigS and yields for the healthy pi = 2.17 
and for the heart failure subject /j,' = 2.4. If we shuffle 
the numbers of the sequence r, we recover 5 = 0.5, a 
fact proving that the memory beyond memory effect is a 
genuine property of heartbeat. 
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ordinal number of heart beat i 
FIG. 2. The inter-beat time T(i) as a function of the 
number of beats, i. The thick line denotes the trajectory cor- 
responding to the coarse graining given by AT = 1/30 sec. 
The vertical lines denote the height of the original data. The 
arrows and the integer labels illustrates how the sequence of 
Tfc's is generated. Inset: Survival probabilities. The circles 
denote the c.h.f. patients and the corresponding fitting func- 
tion is *(i) = 0.19exp[-(i/3.1) 06 ]. The diamonds denote 
the healthy patient and the corresponding fitting function is 
*(i) = 5.71/(0.93 +t) 3 - 25 . 



The additional memory is confirmed by the numerical 
evaluation of the normalized correlation function of the 
variable n — (tj), denoted by C exp (t), where the symbol 
t is the continuous approximation of the discrete patch 
label i. The two- variables model that we are using to 
explain the memory beyond memory effect would yield 



C{t) ex l/t 13 . 



(6) 



The index (3 in that case would be a complicated func- 
tion of the four parameters involved by the two-variable 
model. This means that the DE is a memory detector 
more efficient and much less ambiguous than the corre- 
lation function. The DE selects from the distribution of 
times described by the arbitrary iP(t) the really random 
events yielding a non arbitrary distribution with a unique 



p! . The correlation function C(t), on the contrary, de- 
pends on the details of the model, but does not afford 
an easy way to define them. For the main purpose of 
this letter it is enough to point out that the form of the 
correlation function C exp (t) is 



C exp (t) = {l-e 2 )W(t)+e 2 C(t). 



(7) 



Here W(t) denotes a function dropping from 1 to in one 
time step, while the function C(t), with the asymptotic 
form of Eq. (|6|) , is continuous for for t —* 0. 
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FIG. 3. a) The DE as a function of I. The inset illustrates 
the effect of shuffling, the two lines correspond to 6 — 0.5. 
b) The correlation function C exp (t) as a function of t on a 
healthy and on a c.h.f. individual. 



We account for the structure of Eq.(p as follows. The 
sequence {r^} is generated by the joint use of two mod- 
els. The first is the model of Eq.(|l|) with no additional 
variables and no memory beyond memory property, the 
second is the model with two variables. These two models 
generate two independent sequences {r^}. To any index 
i we assign, with probability e, the value provided by 
the model with additional memory, and, with probabil- 
ity 1 — e, the value provided by the model with only one 
variable. This model is reminiscent of one adopted to ac- 
count for the statistical properties of DNA sequences M . 
The function C exp (t) in one step drops from the value 
Cexp(0) = 1 to the value C exp (\) = e 2 C(l) ~ e 2 , thereby 
allowing us to derive e from the experimental correlation 
function at t = 1. 

In conclusion, the meaning of the parameter e is as fol- 
lows. The value e — 1 would imply that the heart beating 
is described only by a model with two variables, x and y. 
In other words, the larger e the larger the weight of the 



memory beyond memory effect. The parameter 8 is con- 
nected to the time distance between two nearest-neighbor 
real events. If this time distribution is exponential, there 
is no memory of conventional type, as earlier observed. If 
// becomes closer and closer to // = 2, this conventional 
memory becomes stronger and stronger. Thus, to estab- 
lish a more intuitive understanding of what happens to 
memory, regardless of whether it is of the conventional or 
of the new type, let us adopt the following perspective. 
The condition of highest memory corresponds to e = 1 
and 6 = 1. This would mean that the heart beating de- 
pends only on the memory beyond memory model, and, 
that, at the same time, // = 2. The opposite case, of com- 
plete absence of memory, implies e = and // — > oo. This 
would mean that the heart beating is very well modeled 
by the one- variable model of Eq.(gJ), with an exponential 
distribution of waiting times, in other words, without any 
memory of whatsoever form. This leads us to express the 
distribution of patients in the (6, e 2 )-plan of Fig.|| We 
note the surprising result that all the healthy subjects 
and all the heart failure subjects are contained in the top- 
right region and in the bottom-left one, respectively. We 
also notice that all the healthy subjects, but two of them, 
are localized within a small portion of the top-right re- 
gion of the graph, not far from the border with the heart 
failure region. We have the impression that this reflects 
the fact that the healthy function of the heart beating 
system depends on a proper balance of memory and ran- 
domness that the analysis of this letter makes ostensible. 
The distribution of the heart failure subjects within the 
bottom-left region is much broader. It would be desirable 
to have at our disposal the patient survival probability 
as, for example, in Ref. ||, to assess whether a physiolog- 
ical reaction to the c.h.f. pathology, responsible for the 
bottom-left broadened distribution, plays a negative or 
a positive role. The advocates of the second possibility 
might argue that higher randomness and broader distri- 
bution reflect an effort of the perturbed heart beating 
system to explore all possible states to recover the lost 
function. 

It is difficult to establish a connection with the ear- 
lier research work in this field since, to the best of our 
knowledge, the memory beyond memory effect was never 
observed, and we did not find any explicit mention to it 
in the field of heart beating. Is there a correspondence 
between the memory beyond memory effect and the mul- 
tifractal properties observed in Ref. ||? We are inclined 
to believe that there is, due to the similarity between 
Fig.| and Fig.4b of @. 

As to the memory property expressed by 6, and by 
the corresponding //, as well, we would like to mention 
Ashkenazy et al. Ilfl] whose results show for the healthy 
patients a deviation from ordinary scaling higher than 
that of c.h.f patients. If we identify the scaling detected 
by these authors in the case of magnitude fluctuations 
and in the large-time regime, with the scaling parameter 



8 determined by the DE method, we see that our 5 = 
0.82 ± 0.04 for the healthy subjects corresponds to their 
8 = 0.82, and that our 8 = 0.71 ± 0.06 for the heart 
failure subjects corresponds to their 8 = 0.71. Thus, 
we conclude that our findings do not conflict, or do not 
necessarily do, with the earlier findings. 
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FIG. 4. Values of the scaling parameter 8 and of e 2 for 
the healthy and c.h.f. individuals of this analysis. 
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